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Abstract 

The sol-gel transition is studied in two purely entropic models consisting 

of hard spheres in continuous three-dimensional space, with a fraction p of 

nearest neighbor spheres tethered by inextensible bonds. When all the tethers 

are present (p = 1) the two systems have connectivities of simple cubic and 

face-centered cubic lattices. For all p above the percolation threshold p c , the 

elasticity has a cubic symmetry characterized by two distinct shear moduli. 

When p approaches p c , both shear moduli decay as (p — p c )^ where / ~ 2 

for each type of the connectivity. This result is similar to the behavior of 

the conductivity in random resistor networks, and is consistent with many 

experimental studies of gel elasticity. The difference between the shear moduli 

that measures the deviation from isotropy decays as (p — p c ) h , with h ~ 4. 
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Gels are macroscopically large networks, formed when short polymeric units in a solution 
(sol) are randomly cross-linked. The transition from sol to gel is a second order phase 
transition from fluid to solid during which the viscosity of the system diverges and shear 
elasticity (rigidity) develops. Frequently, the geometry of the gels is modeled by percolation 
IjJ, when the monomers are represented by the vertices of some lattice, while the random 
chemical bonds are modeled by bonds joining the vertices with probability p. The gel 
point is identified with the percolation threshold p c , the critical bond concentration above 
which a spanning cluster is formed. Close to p c quantities like the average cluster size 
or the gel fraction have power law dependence on (p — p c ) with exponents independent of 
the details of the lattice. Experimental measurements || of the geometric features of gels 
confirm the correspondence with the percolation model. Therefore, we expect that model 
systems whose geometry is described by percolation will produce a correct description of the 
physical properties of gels, such as elasticity. In this work, we study the elastic properties 
of a purely entropic percolation model of gels, and compare the results with experimental 
measurements of elasticity near the sol-gel transition, as well as with the predictions of 
approximate theories. 

Close to p c the shear modulus of a gel follows a power law: fi ~ (p — p c y . Since the 
polymeric network forming the gel is tenuous and floppy, the dominant contribution to its 
shear modulus is of entropic origin: Upon distortion of the system, the available phase 
space, namely entropy, of the gel decreases, leading to an increase of the free energy, and to 
a restoring force. De Gennes suggested that the exponent / should be equal to the exponent 
t describing the conductivity E of random resistor networks near p c : E ~ (p — p c ) 1 [[J. The 
equality / = t can be proved rigorously for a phantom network [without excluded volume 
(EV) interactions] of Gaussian springs each having the energy E = ^Kr 2 , where r is the 
spring length ||. It describes the entropic elasticity of most phantom networks since the 
latter exhibit an effective Gaussian behavior on sufficiently large scales (provided that they 
are not strongly stretched) ]§,§]. It is an open question whether the equality / = t is also 
valid in the presence of EV interactions. It does, according to some theories which claim 
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that EV interactions primarily influence the bulk (compression) modulus rather than the 
shear (rigidity) modulus of the system, as if the gel is a phantom network embedded into 
a "pressure producing" fluid medium J7j. However, a different approach based on scaling 
arguments concludes that the elastic moduli of a gel are of the order of kT/l; d , where £ is 
the percolation correlation length that diverges as (p — p c ) u , and d is the dimensionality. 
Consequently, the relation / = dp is obtained ||. 

The experimental values of /, measured for different gel systems, are divided between 
the above two approaches. In one group of experiments |§, done on materials like gelatin 
and silica gels, the measured exponent is close to the conductivity exponent t ~ 2 in three 
dimensions [[Hj. Another group consists of experiments in materials like polyester and PVC, 



where the exponent varies from 2.5 to 3.0, and seems to agree with / = dv ~ 2.7 [[Tl|]. The 
gels formed by the materials in both groups of experiments are floppy, and the dominance 
of the entropic contribution to their elastic properties is fairly expected. Thus, the divi- 
sion of experimental works into these two groups is based on the values of the measured 
exponents rather than on the nature of the investigated materials. The origin of the dis- 
crepancy between the experimental results is not clear, and we can only list several possible 
reasons: In some cases the topology of the system does not correspond to three-dimensional 
(3D) percolation model of gels, but is somewhere between gel (cross-linking of monomers 
or short polymeric units) and rubber (cross-linking of a melt of long polymers) [12[ . Addi- 
tional reasons are related to experimental difficulties, such as the imprecise determination 
of concentration of cross-links, or the difficulty to extract the static shear modulus from 
measurements of the low frequency behavior of the dynamic complex modulus. A more 
fundamental reason for the wide range of experimental results is the energetic contribution 
to gel elasticity which mixes with the entropic contribution and influences the "effective" 
exponent. Energetic bending elasticity is characterized by a much larger exponent, / ~ 3.8 
||13|| . Such an exponent is measured only when the entropic contribution to elasticity is 
negligible, e.g, in the experiments in sintered metallic powders [ET3]. When both energetic 



and entropic contributions coexist, we expect the elastic behavior near the gel point to be 



dominated by the latter, since the critical exponent of entropic elasticity (according to both 
approaches to entropic elasticity) is smaller than that of bending elasticity. However, the 
dominance of entropic elasticity near the transition may be limited to a very narrow regime, 
in which the shear modulus is small and difficult to measure. 

In the present work we bypass the problem of mixing of the entropic and elastic contribu- 
tions and investigate the elastic behavior of purely entropic systems: We consider a system 
consisting of hard spheres, connected by "tethers" that have no energy but simply limit the 
distance of a connected pair to be smaller than some value b. We apply a new method in 
which the stress tensor <Jij and the elastic constants CijM of such "hard-spheres-and-tethers" 
systems are measured from the probability densities of contact between spheres and the prob- 
ability densities of having stretched tethers ]15|. The stress and elastic constants are the 



coefficients of the expansion of the free energy density in the components of the Lagrangian 
strain tensor, rjif. f{{r)}) = f({0})+aijr] i j + jC i jkiri i jri k i + . ... In this paper we study systems 
whose elastic properties posses a cubic symmetry under uniform external pressure P. Such 
systems have a diagonal stress tensor: = —P5ij, and only three different non- vanishing 
elastic constants. C\\ = C xxxx = Cyyyy = C zzzz , C%2 = C XX yy = Cyy XX — Cyy ZZ = ■ ■ ■ ] and 



C44 = \{C xyxy + C xyyx ) = . . . [T(|. Cubic systems have two shear moduli /ii = C44 — P and 

// 2 = i(Cii-c 12 )-p. 

The topologies of the networks were defined by considering bond percolation problems on 
simple cubic (SC) and faced-centered cubic (FCC) lattices, with a fraction p of bonds present. 
Each site of the lattice was occupied by a sphere of diameter a, while each present bond 
was replaced by a tether of maximal extension b, which was larger than the nearest-neighbor 
distance bo. Once the topology (connectivity) was defined, the systems were allowed to move 
in a continuous 3D space. For both types of topologies we set the volume fraction of the 
spheres to be 0.2 and the ratio b/a ~ 1.6. We measured the elastic behavior as a function of 
p. It should be noted that even in the absence of tethers a dense system solidifies into FCC 
solid. In such a case we have a first order phase transition, where the liquid-solid coexisting 
phase appears at volume fraction of spheres 0.494, and the system becomes solid at volume 



fraction 0.545 |17| . The topologies of the SC and the FCC systems are quite different: 
In the latter the number of nearest-neighbor lattice sites is larger and, consequently, the 
percolation threshold is smaller: p c ~ 0.12 and p c ~ 0.249 for the FCC and SC topologies, 
respectively. Thus, highly connected rigid regions are formed more rapidly (at lower p) in 
FCC networks. It has been suggested (see, e.g., Devreux et al. in Ref. |9|]) that in real 
gels the creation of such rigid blobs tends to enhance the contribution of energetic bending 
elasticity and, thus, makes the entropy-dominated regime near p c narrower. Therefore, it 
seems interesting to compare the SC and the FCC topologies in a purely entropic model. 
In our Monte Carlo (MC) simulations we used box sizes L = 186 and L = 12V2b for the 
SC and FCC topologies, respectively, with periodic boundary conditions. Figure [l] depicts 
a typical equilibrium configuration of the FCC system. We measured and over a 
broad range of concentrations above p c . Strictly speaking, the rigidity threshold p r is lower 
than the percolation threshold p c due to effects of entanglements |L8j (and, perhaps, also 
due to additional EV effects). However, the two thresholds are so extremely close that they 
are practically indistinguishable in experiments and numerical studies. Therefore, we treat 
p r and p c as identical. The number of quenched topologies and the length of the MC run of 
each individual topology increased as we approached p c . For systems close to p c we needed 
to average the relevant quantities over 10 different topologies, while far from p c , 2-3 sufficed. 
Close to p c the duration of the MC runs is about 500 times larger than the relaxation time 
t of the simulations (see an approximate expression for r in Ref. 0). During each MC run 
the systems were sampled several million times. 

The networks studied in this work posses a cubic symmetry since their topologies are 
defined on cubic lattices. Therefore, their elastic behavior is described by two distinct shear 
moduli rather than one, as in isotropic systems. This property does not exist in experiments 
where the networks are isotropic because of randomness. Figure depicts the two shear 
moduli, /ii and /i2, as a function of {p — p c ) for the SC and FCC systems. For each type 
of connectivity, its own p c is used. The error bars appearing in Figs. |2] and ^ correspond to 
one standard deviation of the averaged quantities. For both systems close to p c , f^i and /i2 



are practically indistinguishable, suggesting that the systems become isotropic. The shear 
moduli can be approximated by the power laws Hi ~ ~ (p — Pc)^ with / = 2.0 ± 0.1 for 
the SC system, and / = 2.1 ±0.1 for the FCC system. Within numerical uncertainty both 
values are similar and consistent with recent estimates of the conductivity exponent in 3D, 
. We already found that / ~ t ~ 1.3 in two dimensions |TIJ and, therefore, we 
expect that / ~ t at any dimension. As we have mentioned earlier in the text, the equality 
/ = t is expected for a phantom network (a = 0) and can be explained by the Gaussian 
nature of its elastic response close to p c . Our results indicate that a similar picture may 
apply to systems with EV interactions. 

In Fig. we observe that the values of [i% and gradually deviate from each other far 
from p c because at large p the systems "remember" the lower (cubic) symmetry of their 
connectivities. For the FCC connectivity fix > //2, while for the SC case \x<i > /J>i- (The 
definitions of the shear moduli /ii and ^ depend on the orientation of the axes of the 
reference system, which in our study were taken along the edges of the conventional cubic 
unit cell.) Figure |3| shows that the difference A/i = |/ix — /x 2 |, follows, in both cases, quite 
similar power laws A/i ~ {p — p c ) h with h = 3.95 ± 0.15 for the SC case and h = 4.15 ± 0.15 
for the FCC case. Because of the similarities of the values of h in SC and FCC systems, 
it is reasonable to assume that h is a new universal critical exponent which characterizes 
deviation from isotropic elastic behavior. While the power law dependence of A/i is not 
surprising due to the self-similar nature of the large percolation clusters, we could only 
support this assumption by numerical data of limited accuracy. We verified the validity of 
the power law dependence on (p — p c ) by attempting (unsuccessfully) to fit the data to other 
functional forms. 

We already saw that the exponent / (describing the leading critical elastic behavior) 
is very similar for self avoiding (SA) and phantom percolating systems. Therefore, it is 
interesting to check whether this similarity applies to the exponent h, as well. For this 
purpose we measured A/i for a phantom FCC bond percolating network with the same 
values of b and bo, but with a = [pPj] . The results of these simulations are also plotted in 
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t ~ 2.0 



Fig. |3|, revealing a power law with h = 4.15 ± 0.15, as in the SA FCC case. The phantom 
Gaussian model which predicts that f — t, cannot be used to predict the value of h since 
it gives Ap = at any bond concentration p P,2l|. Hence, Ap represents deviation from a 
purely Gaussian behavior which originates in the non-Gaussian form of the tether potential 
and (in the SA case) EV interactions. Our results for the exponent h imply that the similarity 
between the critical elasticity of phantom and SA percolating systems may not be restricted 
to the leading Gaussian behavior. 

In conclusion, we have studied the entropic elasticity of 3D purely entropic percolating 
systems. Our study shows that the the critical behaviors of the shear moduli of phantom and 
SA percolation systems are characterized by similar critical exponents which are very close to 
the conductivity exponent. (For phantom systems the elasticity exponent actually coincides 
with the conductivity exponent.) This result agrees with many experimental studies of 
gel elasticity. It corresponds to heuristic theories which assume that the finite clusters 

(1) do not contribute directly to the shear modulus (i.e., behave like a fluid medium) and 

(2) effectively screen out EV interactions in the elastic network. Further support to this 
theoretical description is given by our result for Ap which is also described by similar 
exponents in phantom and SA systems. The exponent h that characterizes the decay of 
Ap seems to be universal, namely independent of the lattice on which the geometry of 
the system is defined, but this point should be established more carefully by studying other 
lattice connectivities, and by measuring h for two dimensional (phantom and SA) percolation 
systems. 
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FIG. 2. Logarithmic plot of the shear moduli \i\ (solid symbols) and ^2 (open symbols) as 
a function of (p — p c ), for FCC (circles) and SC (squares) bond percolating systems. For each 
topology, its percolation threshold is used (p c ~ 0.12 for FCC, and, p c ~ 0.249 for SC). For both 
systems the volume fraction is 0.2 and b/a ~ 1.6. 
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FIG. 3. Logarithmic plot of the difference between the shear moduli Afi as a function of (p— p c ), 
for SA FCC (solid circles), SA SC (squares) and phantom FCC (open circles) bond percolating 
systems. 
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